Monolayer curvature stabilizes nanoscale raft domains in mixed lipid bilayers 
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According to the lipid raft hypothesis, biological lipid membranes are laterally heterogeneous 
and filled with nanoscale ordered "raft" domains, which are believed to play an important role for 
the organization of proteins in membranes. However, the mechanisms stabilizing such small rafts 
are not clear, and even their existence is sometimes questioned. Here we report the observation 
of raft-like structures in a coarse-grained molecular model for multicomponent lipid bilayers. On 
small scales, our membranes demix into a liquid ordered (lo) and a liquid disordered (Id) phase. 
On large scales, phase separation is suppressed and gives way to a microemulsion-type state that 
contains nanometer size lo domains in a Id environment. Furthermore, we introduce a mechanism 
that generates rafts of finite size by a coupling between monolayer curvature and local composition. 
We show that mismatch between the spontaneous curvatures of monolayers in the lo and Id phase 
induces elastic interactions, which reduce the line tension between the lo and Id phase and can 
stabilize raft domains with a characteristic size of the order of a few nanometers. Our findings 
suggest that rafts in multicomponent bilayers might be closely related to the modulated ripple 
phase in one-component bilayers. 



Ever since its introduction some two decades 
ago^l, the lipid raft concept has been discussed 
controversially^, ft re sts on two established facts: (i) 
Biological membranes are laterally heterogeneous. Het- 
erogeneity is necessary to achieve the functions of mem- 
brane proteins, e.g., in cellular signal transduction and 
trafficking or in endocytosis^. (ii) Lipid-lipid phase sep- 
aration is observed in model multicomponent lipid bi- 
layers. A variety of ternary mixtures containing choles- 
terol phase separate into a cholesterol-poor L a or "liquid 
disordered" (Id) phase and a cholesterol-rich "liquid or- 
dered" (lo) phase with a higher degree of chain order®. 
The lipid raft hypothesis states that lipid-lipid phase sep- 
aration contributes to membrane heterogeneity and is 
exploited by nature to organize protein*2H20. Preexist- 
ing "raft domains" supposedly provide a heterogeneous 
environment that sorts proteins and brings them close 
to each other, thus facilitating the protein-protein inter- 
action needed for clustering. This concept provides an 
elegant picture for a number of experimental observa- 
tions, e.g., the reduced mobility of raft-associating pro- 
teins, which depends on cholesterol contentful, the sub- 
microscopic local clustering of raft- associated proteins as 
observed by fluorescence resonance energy transfer^, or 
the existence of detergent resistant membrane fragments 
(DRM) with a high content of cholesterol, sphingolipids, 
and raft-associated proteins^. 

However, all the experimental evidence is rather in- 
direct and the interpretation in terms of lipid rafts has 
been subject to debate^. One problem with lipid rafts 
is that they cannot be observed in vivo with optical mi- 
croscopic techniques, hence they must be tiny. Diffusion 
experiments indicate that the domain sizes of rafts are 
probably in the range of a few tens of nanometers^, and 
they have short lifetimes in the millisecond or microsec- 
ond regime. Thus the current view is that rafts consti- 
tute dynamically changing nanoscale entities, which are 
collected to larger arrays upon need by protein-protein 
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This naturally raises questions regarding the physical 
nature of rafts and the mechanisms stabilizing them. The 
first question is: If the physical basis of rafts is phase 
separation, why are they so small? A number of possi- 
ble explanations have been pointed out. For example, it 
was argued that membranes in vivo are not at thermo- 
dynamic equilibrium, and the constant turnover of lipids 
may well disrupt the formation of large phase separated 
domains®!. Alternatively, it was proposed that immobi- 
lized cytoplasm proteins generate disorder in the mem- 
brane which prevents large scale phase separation 2 ^. 

The second important question is whether biological 
membranes really do tend to phase separate. Model mul- 
ticomponent membranes ma y exh ibit phase separation 
at physiological temperature d 22 * 2 ^ . However, it is not 
clear whether this is also true for biological membranes. 
Veatch et alJ^have isolated giant plasma vesicles directly 
from living rat cells and showed that they do undergo a 
demixing phase transition, but the demixing temperature 
is around T ~ 15 — 25°C. Hence such membranes would 
be in a mixed state at physiological t emp eratures. Veatch 
et al.23 and Honerkamp-Smith et alP3 have argued that 
raft-like structures might emerge as a signature of criti- 
cal fluctuations. This would restrict "rafts" to relatively 
small regions in parameter space, since critical clusters 
only become large close to critical points. 

Alternatively, the membrane might be in the state of a 
two dimensional microemulsion, which is globally homo- 
geneous, but locally phase separated with a characteris- 
tic length scale or domain size^. SchiclP^ has recently 
proposed a mechanism that would stabilize a microemul- 
sion, which builds on a coupling between the local cur- 
vature of the bilayer and the local composition differ- 
ence between the two leaflets 29 . Such a coupling can 
generate modulated phases in mixed membranes under 
tensioipSH^U. Schick argued that it could also stabilize a 
microemulsion-type state with characteristic length scale 
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FIG. 1: Curvature mechanism generating rafts. The liquid 
ordered (lo) and liquid disordered (Id) regions in opposing 
monolayer leaflets are spatially correlated and have different 
spontaneous curvatures. This stabilizes domains of a finite 
size. 



of order 100 nm in membranes under tension, which how- 
ever would diverge in tensionless membranes. 

A more traditional idea is that the membranes con- 
tain line-active agents^ which reduce the line tension 
and eventually turn a phase-separated mixture into a 
microemulsion. Possible candidate agents are proteinic 
or mino rity lipids with a saturated and an unsaturated 
taiOT However, none of the lipids that are considered 
typical for lipid rafts, i.e., cholesterol or sphingolipids, 
have obvious line-active properties^. Thus this mecha- 
nism relies on additional assumptions regarding the com- 
position of membranes containing rafts. 

Typical model membranes for studying rafts ("raft 
mixtures") contain cholesterol and at least two other 
lipid components. Three components seem necessary to 
bring about global lateral phase separation between liq- 
uid membrane phaseiP. However, there is some evidence 
that nanoscopic domains may already be present in bi- 
nary mixtures containing cholesterol, in particular, in 
mixtures of cholesterol and a lipid with high main transi- 
tion temperature (T m ). The literature on these mixtures 
is controversial. Whereas several authors have observed 
immiscible liquid phases, based on various techniques 
such as ESR, NMR, or diffusivity experiments^^, oth- 
ers claim that cholesterol and lipids are miscible in the 
whole high temperature fluid range (see RefP and refer- 
ences therein). In particular, fluorescence microscopy im- 
ages feature only one homogeneous phased FeigensorP^ 
has introduced the notion of "type I" and "type II" mix- 
tures, where the "type II mixtures" exhibit global phase 
separation, whereas "type I mixtures" phase separate 
on the nanoscale, but are globally homogeneous, much 
like microemulsions. The experimental evidence suggests 
that binary lipid-cholesterol mixtures might be type I 
mixtures. Thus bilayers of binary lipid-cholesterol mix- 
tures already seem to have many of the intriguing prop- 
erties attributed to lipid rafts, and a theoretical study 
of such binary systems should provide insight into the 
mechanisms stabilizing rafts. 

In the present paper, we contribute to the raft dis- 
cussion with two main results: First, we present Monte 
Carlo simulations of a coarse-grained molecular model 
for binary lipid bilayers which demonstrate the existence 
of a thermodynamically stable heterogeneous membrane 
phase with raft-like lo nanodomains in a Id environment. 
Hence raft formation is found to be a generic phenomenon 



which can already be observed in binary mixtures and 
does not require specific line-active agents. Second, we 
present a theory which rationalizes our results and ex- 
plains raft formation by a coupling between local compo- 
sition and monolayer curvature. The theoretical picture 
is illustrated in Fig. 1. Liquid ordered domains have a 
propensity to bend inwards. If they oppose each other 

- as suggested by simulations and also by experiments^ 

- the bending competes with bilayer compression. This 
creates elastic tension which reduces the line tension even 
for tensionless membranes, and stabilizes domains with 
a well-defined diameter of the order of 10 nanometers. 

Thus we identify a generic mechanism of raft formation 
in multicomponent membranes. It is different from the 
curvature-mediated mechanism proposed by SchiclP^, 
which is based on bilayer curvature. Both Schick's and 
our mechanism rely on a competition between bending 
and an opposing force. If the bilayer as a whole has 
a propensity to bend, the opposing force is the surface 
tension. Therefore, the surface tension sets the charac- 
teristic length scale which thus diverges in tensionless 
membranes. In our monolayer curvature mechanism, the 
bending is opposed by bilayer compression, which results 
in a characteristic length scale of the order of the mem- 
brane thickness. 



COARSE-GRAINED SIMULATIONS OF 
MIXED LIPID BILAYERS 



Our generic coarse-grained simulation model is based 
on a successful model for one-component lipid bilayers^, 
which reproduces the main phases of phospholipid mem- 
branes (liquid, tilted gel, and ripple phasepOlS and the 
elastic properties of fluid dipalmitoylphosphatidylcholinc 
(DPPC) bilayers at a semi-quantitative leveP^. Here 
we introduce two types of lipids, "phospholipids" P and 
"cholesterol" C , with interactions designed such that C 
is smaller and stiffer than P and has a special affinity to 
P , reflecting the experimental observation that sterols 
in lipid bilayers always tend to be solubilizcd with one 
or two other lipids^. For the purpose of the present 
study, it is essential to design the model such that it 
captures those non-random mixing effects, since they 
most likely drive the local segregation into a Id and a 
lo phase. WaheecP^ and Waheed et alP^ recently stud- 
ied the chemical potential of cholesterol in small systems 
of DPPC/cholesterol bilayers by atomistic (united-atom) 
and coarse-grained simulations. In atomistic simulations, 
they found that the chemical potential drops with the 
cholesterol concentration, indicating a clear tendency of 
local segregation into a cholesterol rich and a cholesterol 
poor phase [about 0.3 thermal enrgy units (kgT) per lipid 
molecule]. This property was not reproduced by coarse- 
grained models that show no sign of random mixing^. 
Our model is described in more detail in the Methods 
section. 

To examine whether the system phase separates lo- 
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FIG. 2: (a) Chemical potential difference in thermal energy 
units (fcsT) vs. cholesterol content from canonical simula- 
tions of small mixed bilayer systems (162 lipids). The corre- 
sponding snapshots show the system in the Id state (b) and 
the lo state (c). The darker chains represent C lipids. 



cally, we first consider small systems (162 lipids) at fixed 
composition. We find that such small systems almost 
always assume one of two states, either disordered (Id) 
or ordered (lo), or jump between the two, depending on 
the composition. Figs. 2 (b) and (c) show sample con- 
figurations of these two states. For future reference we 
have evaluated the pressure profiles across monolayers 
in the two states and computed the spontaneous curva- 
ture of monolayers cq from the first momentP^. In the 
Id state we obtain cq = 0.2 ± 0.2ct -1 , and in the lo state, 
Co = 1.22 ± 0.09rj _1 , using our simulation length unit 
a ~ 0.6 nm. Thus monolayer regions in the lo state have 
a strong tendency to bend inwards, whereas monolayer 
regions in the Id state tend to remain flat. 

The free energy gain fj, for replacing a P chain by a 
C chain (the chemical potential difference) is shown in 
Fig. 2 (a). For temperatures above the main transition 
of the pure P -system (T m = 1.2e), the chemical potential 
curve has an upward slope in the region between w 10% 
and s» 25% C chains. This indicates an unstable regime 
where one would expect spontaneous demixing in larger 
systems. 

When looking at larger systems (20.000 lipids), how- 
ever, we find that the system does not phase separate 
globally. Instead, finite lo domains with a high concen- 
tration of C lipids appear, surrounded by the Id phase al- 
most devoid of C lipids. To ensure that these domains are 
true equilibrium structures and not the result of incom- 
plete phase separation, the simulations were conducted in 
the semi-grandcanonical ensemble at fixed fi, i.e., lipids 
were allowed to switch their identities during the simu- 
lation. Figs. 3 (a) and (b) show an example of an equili- 
brated configuration. One can see the raft-like structure 
of the membrane from the top view (a) and the structure 
of alternating lo and Id regions from the side view (b). 
Consistent with this observation, the behavior of the C 
concentration as a function of [i shows no sign of a phase 
transition [Fig. 3 (c)]. 

We analyze configurations such as shown in Fig. 3 fol- 
lowing an algorithm described in the Methods section. 
This allows us to calculate the distribution of raft sizes 




FIG. 3: Snapshots of a large bilayer system (20.000 lipids) 
featuring raft-like lo domains: (a) top view, and (b) an en- 
larged section of a side view. Parameters are k B T = 1.4, 
fj, a; 6.6 k b T. (c) The chemical potential difference versus C 
( "cholesterol" ) content from semi-grandcanonical simulations 
of large systems. 

(radii of gyration) , shown in Fig. 4 (a) . While the distri- 
butions look rather similar for different values of /i, the 
area fraction (inset) clearly shows that in systems with 
larger [i, i.e., systems with a lower concentration of C , 
a greater fraction of the total raft area is present in the 
form of smaller rafts. 

Next, we address the question whether rafts in oppos- 
ing monolayers are correlated. The normalized cross- 
correlation Ki (see Methods section) between monolayers 
of a configuration i tends to be positive with values rang- 
ing from K — —0.05 to K = 0.15. To analyze whether 
a correlation Ki > or anticorrelation Ki < is sig- 
nificant, we generate a set of configurations with ran- 
domly shifted monolayers and determine the fraction of 
them which have a higher (anti)correlation \Ki\. At large 
fx, i.e., small C concentrations, the distribution of Ki is 
symmetric around zero and shifting monolayers often en- 
hances the correlation. Raft domains are thus uncorre- 
cted in this regime. At large C concentrations, however, 
all values of Ki are positive and shifting monolayers al- 
most always reduces the correlation. We conclude that 
the cross-correlation is significant at higher C concentra- 
tions, and that rafts tend to oppose each other in the 
bilayer. This is compatible with experimental observa- 
tions in membranes exhibiting global phase separation, 
where it was found that lo domains are strongly corre- 
lated across the membrane^. 

Finally, we consider the in-plane structure factor 
[Fig. 5]. For small values of fi we observe a peak at 
nonzero q ~ O.OBa^ 1 , corresponding to a characteristic 
length scale of about 20er ~ 12 nm. Such a peak is a 
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FIG. 4: Characterization of raft domains, (a) Raft size distri- 
bution. The main panel shows the distribution of rafts with 
given radius of gyration. The inset shows the actual fraction 
of the raft area found in rafts of a given size. Lines are guides 
for the eye; ft is given in units of k B T. (b) Cross-correlation 
Ki of configurations i vs. percentile of conformations with 
randomly displaced monolayers, which have a higher correla- 
tion or anticorrelation \Ki\. Every point corresponds to an 
independent simulation configuration i. The more skewed a 
distribution is to the right side, the greater the mean (posi- 
tive) correlation between rafts on both sides. The lower the 
percentile of a point, the less likely this particular value is 
coincidental, a.u., arbitrary units. 
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FIG. 5: Radially averaged structure factor of the raft confor- 
mations for different values of fi, with /x given in units of ksT. 
a.u., arbitrary units. 



typical signature of a microemulsiorP3. 



rafts (see Fig. 1). 

To analyze this possibility, we consider a simple elas- 
tic model for two coupled monolayers with composition- 
dependent spontaneous curvature, which combines a 
model for mixed films by Leibler and AndelmarPSl 
with a bilayer model developed by Dan, Pincus, and 
coworker&SJEU to describe bilayer deformations near in- 
clusions. For simplicity, and in contrast to the model 
proposed by SchidP^, we assume that the local compo- 
sitions on opposing monolayers are strictly equal and do 
not induce bilayer bending. Bending and thickness defor- 
mations then decouple, and for planar membranes, the 
elastic free energy of m onolayer thickness deformations 
can be written as 
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where u(r) denotes the local deviation from the mean 
monolayer thickness to, and the other parameters rep- 
resent material properties of the membrane: the bi- 
layer bending and compressibility modulus k c and k A , 
the spontaneous monolayer curvature c , an associated 
curvature-related parameter Q^, and the Gaussian rigid- 
ity of monolayers ka- Eq. ([I]) holds for tensionless mem- 
branes as well as for membranes under tension^. We 
have used it in the past to fit deformation profiles of one- 
component membranes in the vicinity of inclusions, with 
good results even on molecular length scales235Sl. 

We assume that the monolayer laterally phase sepa- 
rates into two phases, which are separated by narrow in- 
terfaces with a bare line tension Ao- In principle, all mem- 
brane parameters (to, k c , k A , c , C, fcc?) should depend on 
the local composition. For simplicity, however, we will 
assume that only the spontaneous curvature c makes a 
jump from one phase to the other. In that case, the final 
elastic energy after minimization can be written in the 
simple form (see Methods section) 



Fpi = k c Sc n 



dl nVu, 



(2) 



II. ELASTIC THEORY OF RAFT 
STABILIZATION BY CURVATURE 

The main results of the simulations can be summarized 
as follows: (i) lo domains (rafts) of a finite nanoscale 
size are observed in two-component lipid bilayers, (ii) 
they are correlated across the membrane, i.e., rafts on 
the two leaflets tend to oppose each other, and (iii) the 
analysis of pressure profiles in small systems shows that 
lo monolayer domains have a propensity to bend inwards 
(a spontaneous curvature). However, large scale bending 
is prevented by the presence of the rafts on the opposing 
monolayer. These observations motivate the hypothesis 
that the elastic energy associated with the spontaneous 
curvature might be responsible for the finite size of the 



where the line integral J dl runs over all domain bound- 
aries, n is the unit vector normal to the interface, and 
8c a is the curvature mismatch, i.e., the difference of the 
spontaneous curvatures in the inner and the outer phase. 

In many cases, this simple theory can be solved ana- 
lytically. Details of the calculations are presented in the 
Methods section. For isolated plane interfaces, we obtain 
an elastic line energy (an elastic energy per boundary 
length L) 



\% := F e i/L = -$k c Sc 2 /VW^b), 



(3) 



which acts as an additive contribution to the total line 
tension, A t = A + A^ 5 . Note that is negative. 
Here we have introduced the in-plane correlation length 



5 



£ = (fccig/fc^) 1 / 4 and the dimensionless membrane pa- 
rameter b = 2££ 2 /io- Since k c should be roughly pro- 
portional to fc^icpH £ is of the order of the membrane 
thickness. Inserting actual numbers for the fluid phase 
of one-component DPPC bilayers from experiments, all- 
atom simulations, or simulations of our modeP^, one con- 
sistently obtains values around £ ~ (0.9 — 1.4) nm. For 
the membrane parameter &, one obtains b = 0.65 for our 
DPPC model and b = 0.69 for all-atom simulations of 
DPPC. Throughout this paper, we assume |6| < 1. 

Since the elastic contribution to the line tension 
is negative, the effect of elastic relaxation between two 
curvature-mismatched phases is similar to that of adding 
line-active surfactant agents. To assess its impact on the 
demixing transition, we analyze the scaling of Ao and 
X£ close to the critical demixing temperature T c . The 
bare line tension vanishes according to Ao ~ (T c — T) v 
with the critical exponent v = 1, corresponding to the 
universality class of the two dimensional Ising modeP^. 
Likewise, the curvature mismatch Sc will vanish upon 
approaching T c , and it seems reasonable to assume that 
5c is proportional to the composition difference of the 
two phases, i.e., the order parameter of the demixing 
transition. The elastic "line tension" should therefore 
scale as oc Sc 2 Q ~ (T c - T) 2/3 with the 2D Ising ex- 
ponent j3 = 1/8. Comparing the exponents for Ao and 
A^, we find that A^ dominates close to T c , and the line 
tension becomes negative. Thus, macroscopic demixing 
is suppressed at T c . The demixing transition is shifted to 
lower temperatures and gives way to a microemulsion or 
a modulated phase. 

Next, we consider disks with finite diameter D of one 
phase immersed in the other. In that case, the elastic 
line energy depends on D, and we obtain (see Methods 
section) 

el v ) nD 2Vl - b 2 L 2 1 v 2 Ji 

(4) 

where J n and H^P are Bessel and Hankel functions of the 
first kind. For comparison, we also consider the elastic 
line energy for isolated stripe domains of width D. It is 
given by 

Kr(D) = -^=f »(«(! - eMD )- ( 5 ) 

The functional dependence of Xff^(D) and X'^ lp, "(D) on 
the domain size D is illustrated in Fig. 6. As expected, 
both converge towards A^ 3 for large D. At finite D, the 
behavior is nonmonotonic: Both line energies exhibit a 
minimum of similar depth (equal within 3 %), which is 
weak for negative b and becomes more pronounced as b 
approaches b — > 1. The minimum is located at D ~ 4£ 
for disks, and D ~ 2 — 3£ for stripes. This result implies 
that domains of finite diameter become stable even be- 
fore the asymptotic line tension vanishes: Thus the the- 
ory predicts a regime where the membrane is filled with 




FIG. 6: Rescaled elastic contribution to the line tension 
X e i (D), for (a) disk-shaped rafts of diameter D, and (b) stripe- 
shaped rafts of thickness D vs. D, in units of in-plane corre- 
lation length £ and membrane parameter b. 



small nanoscale domains, i.e., rafts. Since disk-shaped 
and stripe-shaped domains have similar line energies, the 
actual shape depends on the composition of the mem- 
brane. At low C concentrations, disks will dominate; 
using £ ~ 1 nm, the predicted characteristic raft size of 
around 4£ corresponds to a few nanometers. This is com- 
patible with typical raft sizes observed in the simulation, 
see Fig. 4. Taking into account that a ~ 0.6 nm, our esti- 
mate is at the lower end of the values typically suggested 
in the literature (~ 10 — 100 nm). 

We should add that the raft domains interact with 
each other, hence the theory actually predicts modulated 
phases with long-range order. Let us consider a system 
where the phase separation - described by a demixing or- 
der parameter $ - is driven by a Ginzburg-Landau free 
energy functional of the form F$ = J d 2 r{|(V<I>) 2 + 
/($)} with /($) = §$ 2 - %<f> 3 + ff$ 4 and let the 
spontaneous curvature Co depend linearly on $ accord- 
ing to c — <&c/£. As shown in the Methods sec- 
tion, a homogeneous phase in this model becomes un- 
stable with respect to modulations with wavelength £ at 
g = 2k c c 2 /(\ — b). Close to the spinodal, we recover 
the Landau-Brazovskii model, which provides a general 
framework for the description of phase transitions driven 
by a short-wavelength instability between a disordered 
phase and ordered phases^!. 

Since the same holds for the models for mixed mem- 
branes with bilayer curvature coupling mentioned in the 
introductio n 1 28 ! 30 ^, we expect the mean-field phase di- 
agrams to be similar, with one important difference: In 
the bilayer coupling case, the characteristic wavelength 
1 / q* \ tends to be in the micrometer range and 

diverges for vanishing membrane tension a, whereas here, 
the characteristic scale £ is in the nanometer range and 
independent of membrane tension. Consequently, the ef- 
fect of fluctuations is expected to be much bigger in the 
present case. Fluctuations are known to shift the order- 
disorder transition and to stabilize a locally structured 
disordered phase via the Brazovskii mechanisrrP 8 . The 
pattern formation then gives way to a microemulsion- 
type raft phase as observed in the simulations. 
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III. DISCUSSION 

To summarize, we have presented a coarse-grained sim- 
ulation model for multicomponent lipid bilayer systems 
containing two types of lipids P and C with properties in- 
spired by phospholipids and cholesterol. Our simulations 
show that this system forms thermodynamically stable 
nanoscale rafts of C -enriched /o-domains surrounded by 
a sea of C -depleted W-phase. The in-planc structure 
factor features a peak at nonzero q, hence we have a 
microemulsion-like structure with a characteristic wave- 
length of around 12 nm. Furthermore, we have suggested 
a mechanism that stabilizes rafts of finite size. The mech- 
anism is based on the idea that the spontaneous curva- 
ture of monolayer regions in the lo- and Id-phase should 
differ. We have established this for our model, and it 
seems likely that it is also the case in real membranes. 
The curvature mismatch then generates elastic interac- 
tions which suppress global phase separation and stabi- 
lize nanoscale structures. 

If rafts are a disordered modulated structure, one 
might ask whether the corresponding ordered modulated 
structure can be observed in nature as well. Indeed, mod- 
ulated ordered structures with a very similar length scale 
are found in one- component membranes of phospholipids 
in the pre-transition region between the fluid phase and 
the tilted gel phase (£3'): The ripple phase PpJ^. It is 
generally observed in lipid bilayers which exhibit a tilted 
gel phase. It is also observed in our modeP^, and the 
periodicity (~ 10 nm) is comparable to the characteris- 
tic length scale of our raft state. The structure of rip- 
pled membranes is much more complicated than that of 
rafts. However, just like in rafts, it involves alternating 
stretches of gel-like and liquid-like domains, and simu- 
lations suggest that ripple formation is to a large ex- 
tent driven by lipid splay, i.e., by monolayer curvature^. 
Thus we may speculate that rafts and ripples represent 
just two sides of the same coin. Curvature-mediated rafts 
might be a generic phenomenon in multicomponent mem- 
branes, just like ripples are a generic phenomenon in one- 
component membranes. 

Unfortunately, rafts are much more difficult to study 
experimentally than ripples due to the lack of long-range 
order, and due to their subsecond lifetimes^. The struc- 
ture factor describing the distribution of domains on the 
nanoscale could possibly be measured by X-ray diffrac- 
tion experiments on aligned multi-lamellar membranes 
("membrane stacks") of the relev ant composition, either 
in reflection or transmissio n 1 61 ! 62 !. This would allow one 
to test whether the in-plane structure factor has a peak 
at nonzero wave vector in type I mixtures as predicted 
by our model [Fig. 5]. In addition, new super resolution 
microscopy techniques^ might provide ways to visualize 
rafts in free membranes on the scale of a few nanometers. 

The curvature mechanism proposed here of course does 
not exclude other mechanisms of raft formation such as 
those discussed in the introduction. Many mechanisms 
might compete in nature. In particular, it should be 



interesting to study the interplay of cu rvature-mediated 
rafts and lipid-protein interaction ;] 21 * 55 ! in future work. 



Methods 

Coarse-grained simulation model 

The model is defined in terms of the length unit a 0.6nm 
and the energy unit e ~ 0.36 ■ 10" 

-20J44] "Phospholipids" (P ) arc- 
represented by simple flexible chains of beads with a hydrophilic 
head and a hydrophobic tail, which self-assemble in the presence 
of structureless solvent beadsEl. "Cholesterol" molecules ( C ) have 
the same basic structure, but they are shorter and stiffer except 
for one flexible end. All lipids are linear chains of six tail beads 
attached to one head bead, connected by finite extension nonlinear 
elastic (FENE) springs with spring constant kf, = 100 ^j, equilib- 
rium bond lengths tq = 0.7<r (P lipid) and r = 0.6cr (C lipid), 
and logarithmic cutoffs at Ar max = 0.2a (P ) and Ar max = 0.15cr 
(C ). Consecutive bonds in a chain with angle © are subject to a 
stiffness potential C/ba(©) = fce(l — cos(B)) with stiffness con- 
stant kg = 4.7e (P lipids), kg = lOOe (C lipids, first four an- 
gles), and kg = 4.7e (C lipid, last angle). Beads that are not 
directly bonded with each other interact via a Lennard- Jones po- 
tential f/ LJ (r/f) = e L j((j:) 12 - 2(i<) 6 ), which is truncated at a 
cutoff radius r c and shifted such that it remains continuous. At 
r c = 1, one recovers the purely repulsive Weeks-Chandler-Anderson 
potentials The interaction parameters for pairs of P or C beads 
(head or tail) and solvent beads are given by 



bead type-bead type 


e/e 




r c /? 


head(any)-head(any) 


1.0 


1.1 


1.0 


head(any)-tail(any) 


1.0 


1.05 


1.0 


head(any)-solvent 


1.0 


1.1 


1.0 


tail(P )-tail(P ) 


1.0 


1.0 


2.0 


tail(P )-tail(C ) 


1.0 


1.0 


2.0 


tail(C )-tail(C ) 


0.9 


1.0 


2.0 


tail(any)-solvent 


1.0 


1.05 


1.0 


solvent-solvent 










Hence all non-bonded interactions except the tail-tail interac- 
tions are repulsive, and the attraction between C tail-beads is 
weaker than that between other tail beads. 

The model was studied by Monte Carlo simulations at constant 
pressure P = 2e/a 3 and constant zero surface tension in a fluctuat- 
ing box of variable size and shape32l. The total number of lipids was 
kept fixed. The composition was sometimes allowed to fluctuate 
(semi-grandcanonical ensemble). In that case, semi-grandcanonical 
moves were implemented by means of configurational bias Monte 
Carlo move during which the identity of a lipid was switched 
between P and C . In the canonical simulations, the same moves 
can be used as virtual moves to determine the chemical potential 
difference u for P and C chains. 



Data analysis 

In order to analyze configurations such as that shown in 
Fig. 3 (a), we map each monolayer onto a discrete grid Xxy, where 
X = 1 stands for "raft" and \ = for "non-raft" . This is done with 
the following algorithm: 

1. Assign chains to upper and lower bilayer leaflet according to 
their head-tail orientation. 
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2. For each layer sort chains into quadratic bins of side length 
1.5(7 according to the icy-Position of the head-bead. 

3. Calculate the number density p, the C density p c and the 
nematic order S for each bin. 

4. At each vertex of the lattice take the mean value of each 
observable in the surrounding bins. 

5. If p ■ p c ■ S > 0.15cr -4 the membrane at the vertex position 
is considered to be in the /o-Phasc. 

6. Apply the Density Based Spatial Clustering of Applications 
with Noise (DBSCAN) cluster detection algorithnffSl with 
Eps=3<T and minimum number of points (MinPts) = 3. The 
clusters are identified as lipid rafts. 

The normalized cross-correlations Ki of opposing monolay- 
ers of a configuration i are calculated according to Ki = 
jf T. x , y (X% ~ X U )(X% ~ X ib )/vito- ib , where x U and X ib are the 
mean values of Xxy an d xLyi respectively, an and are their 
standard deviations, and N is the number of lattice points. To 
analyze whether a correlation Ki is significant, we take the same 
configuration i, displace one of the leaflets by a random offset with 
respect to the other {Xxy = X-x+r m y+r y )< and compare the new 
correlation Ki with Ki. 



Solution of the elastic theory 

The elastic free energy F e i (Eq. jljl) is minimized with respect 
to u(r) in the bulk (inside phase separated domains) and to the 
boundary values of u and the normal derivative nVu at the domain 
boundaries. We note that the latter must be continuous across 
the boundaries, otherwise F e i diverges. This results in the Euler- 
Lagrange equation £ 4 V 4 u + 2b £ 2 V 2 m + u = in the bulk, and in 
boundary conditions at the interfaces between domains: nV(V 2 tt) 
must be continuous, and V 2 « jumps by 2<5c at the boundaries. 
Inserting the Euler-Lagrange equations and the boundary condition 
in Eq. |TJ gives the simplified expression pj. 

For sRipe domains of width D with boundaries at z = ±D/2, 
the deformation profile u(z) satisfying the Euler-Lagrange equation 
and the boundary conditions is given by 



u(z) = 



2C 8c 



,D/2 



s(az)] 



< z < D/2 



VI - b 2 { sR[sin(-aU/2) e iaz ] : z>D/2 



(6) 



with a = \A + *Vl — b 2 /C> and u(—z) = u(z). Two isolated plane 
boundaries are obtained in the limit D — > oo. Finally, the radially 



symmetric analytical solution for disk domains is 
5R[£a Jo(ar) H[ 1} (aD/2)] 



u(r) = 



VI - b 2 



2 



5R[fa MaD/2) Hf'iar)] 



r < R 
r>R 



(7) 



The profile (ml satisfies the boundary conditions by virtue of the 



identity J {z)H\ L '{z) 



.h(z)H^(z) 



-2i/irz. Inserting the ex- 



pressions ^ and ^ into Eq. |2| gives |5j, pj, and (JsJ (in the 
limit D — > oo). 



Derivation of the Landau- Brazovskii model 



To analyze the situation close to the spinodal, we minimize 
F = F e i + F$ for cq = with respect to u. In wavevec- 



tor space q, this gives u q 



4>o 



2££ x((e0 2 ) with X (a:) 



x/(x 2 — 2bx + 1). Inserting this in F and omitting boundary 
terms, we obtain F = | J~J | 'E'q 1 2 <3 2 S'cff (g 2 ) + / d 2 r /($) with 
9off(? 2 ) = g — 4fc c c 2 x((<?£) 2 ) • The function g ff(9 2 ) has a minimum 
at q* =1/5, hence a homogeneous phase (with <3? =constant) be- 
comes unstable at g B n(q*) = 0, i.e., g = g* := 2fc c c 2 /(l — 6). Close 
to the spinodal, contributions q ~ q* dominate, hence we expand 
q 2 g ci f(g 2 ) about q* 2 up to second order. This finally gives a free en- 
ergy expression of the Landau-Brazovskii form F = J d 2 r{ E (A + 



g 2 ) 2 * 2 



2 * 



& + A* 4 } with q% = <7* 2 (1 + (1 - 6)(1 - £)), 



r + \{g - g*)(Qo + 9* 2 ). and T = g* /[2(1 - b)q* 
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